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ABSTRACT 

We investigate the properties of Damped Lyman-a systems (DLAs) using high-resolution 
and large box-size cosmological hydrodynamical simulations of a ACDM model. The nu- 
merical code used is a modification of GADGET-2 with a self consistent implementation of the 
metal enrichment mechanism (Tornatore et al. 2007). We explore the numerical convergence 
of some relevant physical quantities and we vary the parameters describing the properties of 
galactic winds; the initial stellar mass function; the linear dark matter power spectrum and the 
metal enrichment pattern of the IGM (Intergalactic Medium) around DLAs. We focus on the 
properties of dark matter haloes that are likely to be the hosts of DLAs systems: we predict 
relatively low star formation rates 0.01 — 0.1 Af0/year) and metallicities around 0.1 Zq, 
at least for the bulk of our haloes of masses between 10^ and 10^*^ H^^Mq hosting DLAs. For 
more massive haloes metallicities and star formation rates depend on the specific wind model. 
We find that strong galactic winds with speed of about 600 km/s, in an energy-driven wind 
scenario, are needed in order to match the observed column density distribution function for 
DLAs and the evolution of the neutral hydrogen content with redshift. The momentum-driven 
implementation of the galactic wind model, that relates the speed and mass load in the wind 
to the properties of the dark matter haloes, shows a behaviour which is intermediate between 
the energy-driven galactic winds of small 100 km/s) and large 600 km/s) velocities. 
At z = 3 the contribution of haloes of masses between 10^ and 10^° H^-^Mq, for DLAs 
below 10^"'^ cm^-^, to the column density distribution function, is significant. By interpo- 
lating physical quantities along line-of-sights through massive haloes we qualitatively show 
how different galactic wind models impact on the IGM around DLAs. Furthermore, we anal- 
yse statistics related to the velocity widths of Sill associated to DLAs: while the expanding 
shells of gaseous matter associated to the wind can account for the observed velocities, the 
metallicity in the wind seems to be rather clumpy and this produces an underestimation of the 
observed velocity widths. We outline possible solutions to this problem. 

Key words: cosmology: theory - intergalactic medium - galaxies: formation - quasars: ab- 
sorption lines - methods: numerical 



1 INTRODUCTION 

Damped Lyman-a systems (DLAs) are defined as quasar (QSO) 
absorption syste ms with neutral hyd rogen column density A'^hi > 

2 X 10^° cm"^ ( |Wolfe et aLl|l986| . DLAs are considered as an 
important reservoir and/or sink of gas for the galaxy formation pro- 
cess in the high redshift universe and their neutral hydrogen content 
dominate the total neutral hydrogen budget over a large fraction of 
the cosmic history. The interplay between DLAs and galaxies is 
thereby fundamental and should be addressed by any galaxy for- 
mation model. 

Significant observational efforts have been made in order to 
understand the nature of DLAs (e.g. |Wolfe et aL]|1995[ |Storrie-| 



ILombardi & Wolfe|2000l |Rao & Turnshek|2000| [Prochaska et al.| 
"ZSOT; "Per oux et all [20031 |Chen & Lanzetta] [2003] |Prochaska| 
et al. 2005p , while, on the theoretical side, semi-analytical mod- 
els and high-resolution hydrodynamical simulations are routinely 
performed to investigate the relation between DLAs and dark mat- 
ter haloes and to match their observed properties (e.g. |Katz et al 



p96; Gardner et al. 1997, Haehnelt et al. 1998, Maile r et al. 



12001; Gardner et al. 2001, Nagamine et al. 2004, Okoshi & Na- 



gashima|2005[ [Nagamine et al.|2007||Pontzen et al.|2008[ [Barnes 
|& Haehnelt|2008^ . Despite these efforts, the nature of DLAs is still 
unclear (for a review see [Wolfe et al.|2005| l and the interpretation of 
the large observed velocity width of low-ionization species is not 
unambiguous. In fact, the disk-like model of |Prochaska & Wolfe[ 
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([T997J in which DLAs are thick rotating disks with speed typical 
of a present-day spiral galaxies and the alternative model based 
on the assumption that DLAs are protogalactic clumps (Haehnelt' 
|et al.|1998 l seem to be both viable. More precisely, Haehnelt et al.| 
1 1998^ showed that while the velocity width profiles of DLAS can 



be reproduced by rotation in disks as well as protogalactic clumps, 
the velocity width distribution cannot be reproduced by rotation in 
disks in the context of a cold dark matter model. 

In order to make significant progress in this field observational 
and theoretical/numerical efforts are needed to tackle the problem 
under many different aspects. Recently, the SDSS (Sloan Digital 
Sky Survey) has allowed to measure to an unprecedented precision 
the statistical properties of DLAs (incidence rates, column density 
distribution function, etc.) over a wide redshift range ( |Prochaska| 
|et al.|2005] >. Moreover, attempts to find the galaxy counterparts of 
DLAs represent the most promising way to understand their nature 
(e.g. |Fynbo et al.|7999 ; Christensen et al. 200 71. High-resolution 
spectroscopic studies have also played a crucial role since low- 
ionization metal lines, that are supposed to be good tracers of the 
neutral gas in DLAs, can be identified shedding Ught on the DLAs 
chemical and physical properties (Matteucci et al."1997; Calura 
et al. 12003 1 [Vladilo & Pgroux 2005 ; Vladilo et al. 2006 . Prochaska 



& Wo lfe 2008] [Vladilo et al.|2008^ . Here, we choose to focus on 



the role of feedback in the form of galactic outflows and its impact 
on the DLAs properties. 

Recently, Pont zen et al.| ( |2008] l analysed several different hy- 
drodynamical simulations with and without feedback using an ap- 
proximate a-posteriori radiative transfer scheme. They matched 
most of the observed DLAs properties apart from a tension with 
the observed velocity widths. In their framework, haloes of virial 
masses between lO" and lO" h'^MQ were the main contributors 
to the DLA cross-section. [Barnes & Haehnelt] p008| l with a semi- 
analytical model claimed instead that in order to reproduce the ve- 
locity widths distribution the cross-section of haloes less massive 
than 10^" Mq should be exponentially suppressed. Observa- 
tions of the HI distribution in the local universe made by |Zwaan| 
|et al.|p008| l show that the link with galactic superwinds might be 
stronger than expected and favoured this explanation instead of the 
disk-like model to explain the observed velocity widths. Also some 
correlation properties of neutral hydrogen rich galaxies and HI ab- 
sorbers in the local universe seem to be better fit when galactic 
winds are taken into account d Pierleoni et al.|2008^ . However, the 
situation is quite confusing and potential problems for this interpre- 
tation are extensively discussed in |Prochaska et al.]^2008^ . 

In this work, we extend the analysis performed recently by 
|Nagamine et al.|(2007[ l that addressed the properties of DLAs in a 
ACDM universe. The main differences can be summarized as fol- 
lows and will be more extensively discussed later: i) we rely on 
simulations that have on average a factor 10 better mass resolution; 
ii) we use a different version of the hydrodynamical code GADGET- 
2 that incorporates the dependence of the radiative cooling function 
on the global metallicity of the gas (following [Sutherland & Do^ 
pita|1993] l and a self consistent metal enrichment model; iii) we 



explore the effect of varying the stellar initial mass function and 
most of the parameters describing the wind model. As for the wind 
model we focus mainly on energy-driven scenarios but we allow for 
an extra simulation that is based on a momentum-driven wind. At 
least for the low-density IGM, this model seems to better fit some 
observational properties ( e.g. [Oppenheimer & Dave|2006[[2008] l. 
The goal is to see which DLAs properties can be reproduced by the 
hydrodynamical simulations and to investigate closely the impact 
of galactic winds both on the neutral hydrogen and on the metal 



distribution around galactic environments in the high redshift uni- 
verse (mainly in the range 2 = 2 — 4). 

This work is intended as a preliminar quantitative attempt to 
match some observed properties of the metal and neutral hydrogen 
distribution using a self consistent chemo-dynamical code that has 
been already tested for the intra-cluster medium ( [Tornatore et al.[ 
[2007| >. After having addressed some global properties and evolution 
of IGM low-ionization species, we will focus on DLAs because 
their statistical properties are well measured over a wide redshift 
range and they could provide a useful benchmark for the physics 
implemented in our simulations. 

This paper is organized as follows. In Section 2 we describe 
our set of simulations along with the two different galactic feed- 
back implementations used: energy-driven and momentum-driven 
winds. In Section 3 we compare some global properties of the simu- 
lations: in particular temperature and metallicity relations (Section 
3.1), star formation rates and evolution of ion species (Section 3.2). 
In Section 4 we focus on the neutral hydrogen distribution around 
galaxy-sized haloes and we study DLAs properties like the cross- 
section (Section 4.3), the incidence rate (Section 4.4) and the col- 
umn density distribution function (Section 4.5). Section 5 is dedi- 
cated to the simulated QSO spectra extracted from our simulations, 
to study the velocity width distribution of low-ionization species 
(Section 5.1) and the correlation between metallicity and velocity 
widths (Section 5.2). Finally in Section 6 we summarize our main 
results and we draw some conclusions. 



2 THE SIMULATIONS 

We run a number of hydrodynamical simulations in order to ex- 
plore the effect of changing box-size, numerical resolution, stellar 
initial mass function, prescription for energy feedback and nature 
of dark matter. We use a modified version of GADGET-2 , a parallel 
Tree-PM Smoothed Particle Hydrodynamics (SPH) code ( Springe!] 
[2005[ >. The main modification consists in an accurate modelling of 
the chemical evolution which allows us to follow the metal release 
from Type II and Type la supernovae (SNII, SNIa), along with low 
and intermediate mass stars (LIMS). We refer to [Tornatore et al.] 
( [2007) ; [Fahg"an et al.[ ( |2008) ; [Borgani etaTlpOOSl l for more details 
while we briefly outline in the following the most important fea- 
tures of the chemical enrichment model. 

The simulations cover a cosmological volume (with periodic 
boundary conditions) filled with an equal number of dark matter 
and gas particles. The cosmological model chosen is a flat ACDM 
with the following parameters: Qom ~ 0.24, f^ob = 0.0413, 
Qa = 0.76, Us = 0.96, Ho = 73 km/s/Mpc and ag = 0.8, 
which are in agreement with the latest results from large scale struc- 
ture observables such as the cosmic microwave background, weak 
lensing, the Lyman-a forest and the evolution of mass function 
of galaxy clusters ( Komatsu et al.|2008{ [Lesgourgues et al.|2007[ 



[Vikhlinin et al.[2008| . The input linear dark matter power spectrum 
for the initial conditions has been generated at z = 99 with CMB- 
FAST ( Sel jak & Zaldarriaga|1996"^ and includes baryonic acoustic 
oscillations. In one case only we change the initial linear dark mat- 
ter power spectrum by running a warm dark matter (WDM) simu- 
lation with a modification (suppression at the small scales) of the 
initial power spectrum to account for a thermal dark matter particle 
of mass of 1.2 keV, which is in rough agreement with recent results 
obtained from Lyman-a high-resolution QSO spectra at redshifts 
similar to those investigated here ( [Viel et al.|2008[ l. 

Radiative cooling and heating processes are included foUow- 
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ing the implementation of |Katz et al.H1996) . We assume a mean 
Ultra Violet Background (UVB) produced by quasars and galaxies 
as given by |Haardt & Madau| ( |1996[ l, with the heating rates multi- 
plied by a factor 3.3 in order to better fit observational constraints 
on the temperature evolution of the Intergalactic Medium (IGM) at 
high redshift. This background gives naturally ar_i2~0.8 — lat 
low redshift in agreement with recent observational measurements 
( [Bolton et al.|2005|[Faucher-Giguere et al.pOOS ^. Multiplying the 
heating rates by the factor above (chosen empirically) results in 
a larger IGM temperature at the mean density which cannot be 
reached by the standard hydrodynamical code but aims at mimick- 
ing, at least in a phenomenological way, the non-equilibrium ion- 
ization effects around reionization (see for example [Ricotti et al.| 
[2000)|Schaye et al.|2000||Bolton & Haehnelt|2007| >. We follow self 
consistently the evolution of the following elements: H, He, C, O, 
Mg, S, Si and Fe. The contribution of metals is included in the cool- 
ing function adopting the tables of |Sutherland & Dopita| ( |1993] l, that 
assume the solar value for the relative abundances. In this paper we 
use the solar metallicity and element abundances given by |Asplundl 
[etaLl ( |2005t . 

The standard multiphase star formation criterion is used in 
which an effective prescription for the inter-stellar medium (ISM) 
is implemented jSpringel & Hernq uist 2003). In this effective 
model, a gas particle is flagged as star forming whenever its den- 
sity exceeds a given density threshold, above which that particle 
is treated as multiphase. With this prescription baryons are in the 
form either a hot or a cold phase or in stars, thereby this density 
threshold marks the onset of cold clouds formation. Following |Tor-| 



natore et al. 



(2007 



we set the threshold value to pth = 0.1 cm" 
in terms of the number density of hydrogen atoms. 

The neutral hydrogen fraction /hi is associated to each gas 
particle and is stored in each simulation snapshot. However, we 
follow [Nagamine et al.| ( |2004) l to assign a-posteriori a new mass in 
neutral hydrogen to gas particles above the density threshold which 
reads: 

rriYii = fmrnn {p < pth) (1) 

mm = fcmn (p ^ pth) , (2) 

with /hi the neutral hydrogen fraction that depends on the UVB 
used, rriH the hydrogen mass of the particle (/hi and rriu are deter- 
mined self consistently inside the code), /c the fraction of mass in 
cold clouds and pth the star formation threshold. Here fc = Pc/p, 
where pc is the density of cold clouds and p the total (hot -I- cold) 
gas density. Individual molecular clouds cannot be resolved at the 
resolution reachable in cosmological simulations, thus pc repre- 
sents an average value computed over small regions of the ISM. 
We refer to |Springel & Hernquist|pb03^ to better understand how 
this cold fraction is related to the physics of the ISM. 

Besides including different contributions from SNII, SNIa and 
LIMS, our model of chemical evolution accounts for the age of 
various stellar populations. Metals are thereby released over differ- 
ent time-scales by stars of different mass. For the stellar yields we 
use: SNIa -|Thielemann et al.|(j2003); SNII - Woosley & Weaver] 
( [1995) ; LIMS - van den Hoek & Groenewegen ( 1997). The mass- 
range for the SNII is m > 8A/q, while for the SNIa is m < 8Mq 
with a binary fraction of 10%. We also adopt the lifetime func- 
tion given by 'Padovani & Matteucci| ( |1993) . Finally we use three 
distinct stellar inital mass functions (IMFs): a Salpeter, a Kroupa 
and an Arimoto-Y oshii IMF . Our re ference choice is the functional 
form proposed by Salpeter jl955| : ^{m) oc m~^, where ip{m) 
is the IMF by mass and y = 1.35. Arimoto-Yoshii ( [Arimoto~& 
|Yoshii|[T987l ) IMF has ^ay = 0.95, while Kroupa IMF jKroupa 



|2001[ l adopts a multi-slope approximation: j/kh ~ 0.3 for stellar 
mass m < O.5M0, j/kr = 1.2 for O.5M0 5C m < IMq and 
J/KH, = 1.7 for m ^ IMq. 

2.1 Energy-driven winds 

The implementation of galactic winds used is extensively described 
in Springel & Hernquist ( 2003 1 however, since this paper addresses 
the effects on the chemical enrichment of winds blowing from high 
density regions, we summarize here the major details of the original 
model and the modifications that we made. 

In the multiphase model of |Springel & Hemquist|p003| l gas 
contained in dark matter haloes cool and settle into rotationally sup- 
ported discs where the baryons are converted into stars. However, 
in the star-forming multiphase medium, it is plausible that not all of 
the hot gas in supernova remnants remains confined to the disc: su- 
pernova bubbles close to the surface of a star-forming region may 
break out of a disc and a galactic-scale wind associated with star 
formation may develop. 

Basically, the wind mass-loss rate M„ is assumed to be pro- 
portional to the star formation rate M* according to J\fw ~ 
where 77 is the wind efficiency. Star-forming gas particles are then 
stochastically selected according to their star formation rate to be- 
come part of a blowing wind. Whenever a particle is uploaded to 
the wind, it is decoupled from the hydrodynamics for a given period 
of time (calculated as mentioned below) or till the density around 
drops below a given fraction of the density threshold set for the 
onset of the star formation, in order to effectively reach less dense 
regions. In our case, this limiting density for decoupling is 0.5pth. 
This allows the wind particle to travel 'freely' up to few kpc until 
it has left the dense star-forming phase, without directly affecting 
it: only outside the disc the hydrodynamical interactions within the 
halo could stop the wind. Unlike in [Springel & Hemquist|p003| ), 
we decide here to fix the velocity of the winds, v^, instead of the 
fraction of the energy made available by SNII explosions to power 
galactic ejecta. Thus, four parameters fully specify the wind model: 
the wind efficiency rj, the wind speed Hw, the wind free travel length 
Zw and the wind free travel density factor (5w. 

The maximum allowed time for a wind particle to stay hydro- 
dynamically decoupled is then t^cc ~ Ivi/vvi- The parameter 
has been introduced in order to prevent a gas particle from getting 
trapped into the potential well of the virialized halo and in order 
to effectively escape from the ISM, reach the low density IGM and 
pollute it with metals. [Nagamine et al.| ( (2007^ showed that global 
DLAs properties are relatively insensitive to the value of l^. 

In this work, we will consider two values for the wind velocity: 
Ww = 100 km/s (weak winds, WW) and u„ = 600 km/s (strong 
winds, SW). In our implementation the parameter rj is kept fixed 
to the value 2. For a more quantitative description of the numerical 
implementations made and for an application to the intra-cluster 
medium physics we refer to the paper by |Tomatore et al.[j2007^ . 

2.2 Momentum-driven Winds 

In this Section we describe our implementation of the momentum- 
driven galactic winds, which have been first implemented in simu- 
lations by Oppenheimer & Dave|f2006||200 8). These authors sug- 
gested that with this scheme of winds simulations reproduce the 
statistics of CIV absorption in the high-redshift IGM better than 
using the above energy-driven winds. 

Observations of outflows from starbust galaxies have im- 
proved considerably in recent years. |Martin| ( |20()5) found that the 
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Run 


Size (h-^Mpc) 


A'gas 




soft, (h ^kpc) 


Wind (km/s) 


IMF 


Zf 


SW 


10 


320^ 


3.5 X IQS 


1.5 


600 


Salpeter 


2.25 


WW 


10 


320=* 


3.5 X IQS 


1.5 


100 


Salpeter 


2.25 


MDW 


10 


320=* 


3.5 X 10^ 


1.5 


cr-dependentf"' 


Salpeter 


2.25 




10 


320^ 


3.5 X 10^ 


1.5 


600 


Ki'oupa 


2.25 


SWay 


10 


320^ 


3.5 X IQS 


1.5 


600 


Arimoto-Yoshii 


2.25 


SWio,448 


10 


448^ 


1.2 X IQS 


1 


600 


Salpeter 


3.00 


SWl0,256 


10 


256^ 


6.8 X IQS 


2 


600 


Salpeter 


2.25 


SW20,512 


20 


5123 


6.8 X IQS 


2 


600 


Salpeter 


3.00 


SW5,320 


5 


320^ 


4.3 X 10" 


0.75 


600 


Salpeter 


3.00 


SWwDM^*") 


7.5 


320^ 


1.5 X 10^ 


0.75 


600 


Salpeter 


3.00 



Table 1. Summary of the different runs. Column 1, run name: SW, Strong Winds, WW, Weak Winds and MDW, Momentum Driven Winds; column 2, 
comoving box size; column 3, number of gas particles; column 4, mass of gas particle; column 5, Plummer-equivalent comoving gravitational softening; 
column 6, wind speed; column 7, Initial Mass Function (IMF) chosen (see text); column 8, redshift at which the simulation was stopped, (a): a is the velocity 
dispersion of the halo that host the "wind" particle (see Section [2!2) . (6): WDM (Warm Dark Matter) simulation with a modification (suppression at the small 
scales) of the initial linear dark matter power spectrum to account for a thermal dark matter particle of mass 1 .2 keV. 



terminal wind velocity scales roughly linearly with circular veloc- 
ity, with top winds speeds around three times the galaxy's circular 
velocity. [Rupke et aL| ( |2005l l studied a large sample of luminous in- 
frared galaxies and found that, at least when combined with smaller 
systems from |Martin| ( |2005) , those trends continue to quite large 
systems. 

A feasible physical scenario for the wind driving mechanism 
is derived by noting that the observed scaling are well explained 
by a momentum-driven wind model such as that outlined by Mur- 
|ray et aT] \2QQ5) . In such a scenario, the radiation pressure of 
the starbust drives an outflow, possibly by transferring momen- 
tum to an absorptive component (such as dust) that is hydrody- 
namically coupled with the gas component, which is then dragged 
out of the galaxy. Following Oppenheimer & Dave ( 2006 2008 1 
we test a single model (very close to their "mzw" run) based on 
momentum-driven winds. In such a model the wind speed scales 
as the galaxy velocity dispersion a, as observed by I Martin] p005 | l . 
Since in momentum-driven winds the amount of input momentum 
per unit star formation is constant, this implies that the mass load- 
ing factor must be inversely proportional to the velocity dispersion. 
We therefore use the following relations: 



vind ~ 3cr \/7l 



(3) 



where /l is the luminosity factor in units of the galactic Edding- 
ton luminosity (i.e. the critical luminosity necessary to expel gas 
from the galaxy), cro is the normalization of the mass loading fac- 
tor and we add an extra 2a kick to get out of the potential of 
the galaxy in order to simulate continuous pumping of gas until 
it is well into the galactic halo. The functional form of the lumi- 
nosity factor includes a metallicity dependence for /l owing to 
more UV photons output by lower-metallicity stellar populations: 

/l = /l,0 X lQ"0"029(l°gZ+9)"^+0-417694^^1^g^g^^^jij^,2005l 

suggests: /l,© « 2. The mass loading factor controls star forma- 
tion at early times, so cto can also be set by requiring a match to 
the observed global star formation rate. Following |Oppenheimer &| 
[Dave (2008) we set cro = 150 km/s. We determine cr directly from 
the siinulation only for the most massive haloes, while for haloes 
with masses below 1.7 x 10^ h~^MQ (i.e. formed by less than 10'^ 



dark matter particles) we use the relation calibrated by |Evrard et al.| 
( |2008^ , using a variety of N-body simulations: 



0"Dm(A^, 2) ~ CrDM,15 



h{z)M2oo 



IOIsMq 



km/s. 



(4) 



where ctdm.is = 1082.9 ± 4.0 km/s is the normalization at 
mass lO^^/i"^M0, a = 0.3361 ± 0.0026 is the logarithmic 
slope and M200 is defined as the total mass within a sphere with 
mean interior density 3M200 /47rr2oo = 200pc{z) {pc{z) is the 
critical density at redshift z). Even for this wind implementation 
the particles are stochastically selected in the same way as for the 
energy-driven scenario. 

In Table I we summarize the main parameters of the cosmo- 
logical simulations performed including the mass associated to the 
gas particles and the gravitational softening. All the simulations 
start at redshift z = 99. The reference simulations are SW ("Strong 
Winds"), WW ("Weak Winds") and MDW ("Momentum-Driven 
Winds") with a total of 2 x 320^ dark matter and gas particles 
in a 10 comoving h^^Mpc box and Salpeter IMF. Furthermore, 
SWkr (Kroupa IMF) and SWay (Arimoto-Yoshii IMF) explore 
the effects of a different IMF compared to the Salpeter one. Some 
simulations are intended to explore resolution and box-size effects 
like SWio,448, SWio,256, SW2o,5i2 and SW5,320. The most CPU- 
time consuming run is the SW2o,5i2, which significantly increase 
the statistics of dark matter haloes able to host DLAs. Finally, the 
SWwDM aims at exploring the effect of a modification in the initial 
linear dark matter power spectrum. 



3 GLOBAL PROPERTIES OF THE SIMULATIONS 

3.1 IGM temperature and metallicity relations at 2: = 3 

In this Section we investigate two global properties of the IGM, 
namely its thermal and chemical state, focussing on the differences 
between the SW, WW and MDW runs at 2 — 3. A more detailed 
analysis of global IGM properties down to small redshifts is beyond 
the scope of this paper and will be presented in a future work. 

In Figures [T] and |2] we show the IGM temperature-density re- 
lation and the metallicity-density relation for the WW (left panel), 
SW (middle panel) and MDW (right panel) simulations. In these 
figures it is shown an "effective" temperature for the star-forming 
particles determined as the mean temperature weighted by the con- 
tribution in mass of the cold and the hot phases (see Section 2). In 
all cases the phase diagrams are color coded according to the gas 
mass fraction, as indicated by the vertical bars. For all the three runs 
a significant fraction in mass (about 10%) of the IGM resides in the 
tight power-law relation T = To(l + 5)''~^ at around the mean 
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-2 02458 -2 02468 -2 02468 

Log p/<p> Log p / <p> Log p / <p> 



Figure 1. Gas (IGM) temperature-density relation (without explicitly splitting the star-forming pailicles in the hot and the cold phases), for the WW (Weak 
galactic energy-driven Winds of 100 km/s; left panel), SW (Strong galactic energy-driven Winds of 600 km/s; middle panel) and MDW (Momentum-Driven 
galactic Winds; right panel) simulations at z = 3. The vertical bar indicates the gas mass fraction. 




Log p/<p> Log p / <p> Log p / <p> 

Figure 2. Gas (IGM) metallicity-density relation (without explicitly splitting the star-forming particles in the hot and the cold phases), for the WW (Weak 
galactic Winds of 100 km/s; left panel), SW (Strong galactic Winds of 600 km/s; middle panel) and MDW (Momentum-Driven galactic Winds; right panel) 
simulations at z = 3. The vertical bar indicates the gas mass fraction. 



density. Gas at this density is usually responsible for Lyman-a for- 
est absorption. At 2 = 3 the IGM has a temperature at the mean 
density of about To = 15000 K and the slope 7 of the temperature- 
density relation is around 1.6. In our simulations we assume pho- 
toionization equilibrium and under this assumption it is difficult 
to get values of 7 < 1.6, which however are in better agreement 
with the rece nt measurements of the Lyman-Q probability distribu- 
tion function ( [Bolton et al.|2008| l. An IGM fraction between 10"" 



and 10~ is in the form of shock-heated gas at temperatures higher 
than 10^ K. The tail of cold gas at large densities resides in the 
inner parts of the haloes while at even larger densities, > 10*(p), 
and temperature, > 10* K, is associated to star-forming regions. 
We note a plume of cold and overdense gas particles that are below 
10* K and carry a fraction of mass of the order 10"'': these are gas 
particles that are in the wind phase. After these particles enter the 
wind regime they become hydrodynamically decoupled for a cer- 
tain period of time (see Section 2.1). During this period they travel 
'freely' towards regions in which both the density and the tempera- 
ture are lower. As a consequence, when these particles re-couple to 
the hydrodynamics, their temperature drop because of the adiabatic 
expansion. 

In Figure [2] we show the metallicity-density relation for the 
three simulations. In the SW case the gas that has been enriched 



spans a wide range of densities and can attain values that are very 
close to solar. Also gas that is below the mean density appears to 
be metal enriched at a high level in this simulation. The bulk (10"* 
in mass fraction) of the IGM is either in the form of gas at the 
mean density or in the form of very dense gas at about 10* (p) with 
metallicities of about 0.1 Zq. 

All the previous findings are particularly interesting when the 
WW (left panel) and SW (middle panel) runs are compared. For 
the WW temperature-density relation it is clear that the amount of 
shock-heated gas is signifcantly reduced and the region of the di- 
agram at temperatures between 10^ and 10^ K is less populated. 
Thereby, strong galactic winds with speed of order 600 km/s heat 
the IGM significantly, while weaker winds with speed of order 100 
km/s are less effective in doing this. The gas particles in the wind 
phase at T < \Q*K are now fewer than in the SW case and at 
higher overdensity. Moreover, the star-forming high density tail is 
considerably more extended than in SW case confirming the less 
efficient feedback of WW in suppressing star formation. In con- 
trast with the SW simulation, in the WW case there is slightly less 
metal-enriched underdense gas, while the metallicity gradient with 
density is steeper. Higher values of metallicities are now reached, 
which attain super-solar values at very high density. In the WW 
case the bulk of the enriched gas shifts from the mean density val- 
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Figure 3. Left Panel: cosmic star formation rate (SFR) for some of the hydrodynamical simulations of Table I. Right Panel: evolution of the total Cjjl (x 1000) 
as a function of redshift for some of the hydrodynamical simulations of Table I. 



ues of the SW case to a region which is about 10* times denser than 
the mean. 

We stress that the difference between the SW and WW simu- 
lation is only in the wind speed which is 600 and 100 km/s, respec- 
tively, while all the other parameters have been kept fixed. These 
values of the wind velocity have been chosen in order to embrace 
values suggested by observational studies of HI and CIV-galaxy 
correlation and Lyman-a forest in close QSO pairs ( [Adelberger| 
|et aI.|2003]|Rauch et al.|2005| l, and observations of interstellar lines 
in ultraluminous infrared galaxies (Martin 2005 ). 

The right panels of Figures[T|and|2]show results for the MDW 
simulation. Looking at the temperature-density diagram one can 
easily see a general trend that is intermediate between those of 
the SW and WW runs. For example, the region associated to wind 
particles (cold and overdense particles below 10'* K) is less pro- 
nounced with respect to SW but not as negligible as the WW run. 
However the metallicity-density relation is quite different from the 
energy-driven implementation. In particular the metallicity reaches 
values slightly higher than SW and WW for the whole range of den- 
sities and most importantly there is much more enriched material 
around the mean density. This result suggests that MDW is more 
efficient in polluting the IGM above and around the mean density: 
this is due to the fact that, unlike for SW and WW, in MDW even 
small haloes contribute significantly to the enrichment. In fact, such 
small haloes have small velocity dispersion and correspondingly 
greater loading factors (see Eq. (|3|) for their winds. 



3.2 Star formation rates and evolution of ion species 

In this Section we analyse the star formation rates (SFR) for the dif- 
ferent simulations and the evolution of the neutral hydrogen content 
and of two ion species, CIV and OVI, that are usually observed in 
absorption in QSO spectra even at higher redshifts than those ex- 
plored here (e.g. jBecker et al.|2008[ [Ryan- Weber et al.|2009^ . 

In the left panel of Figure [3] we show the total star forma- 
tion rate of the simulated volume as a function of redshift. The 



WW simulation has nearly a factor 10 higher star formation rate 
compared to the others at 2 ~ 3. This means that the feedback 
mechanism induced by strong galactic winds is effective in pau- 
pering the metal-rich star-forming gas significantly decreasing the 
star formation rate. The SW2o,5i2 and SW5,32o embrace all the SW 
runs and this is due to cosmic variance effects. We note negligible 
differences between the Salpeter IMF and the other two IMFs used 
(Kroupa and Arimoto-Yoshii), at least at the relatively high red- 
shifts considered here. We underline, however, that we neglect the 
effect of assuming different IMFs on the observationally inferred 
cosmic star formation rate; this means that we do not change the 
star formation efficiency as it would be required in order to match 
the observables when the number of massive stars per unit mass 
of formed stars changes. That is because here we are more inter- 
ested in the chemical and energetic effect of the IMF. We also point 
out the intermediate trend of the MDW, in between SW and WW. 
At very high redshift {z > 10) the SW5,32o and SWio,448 show a 
higher star formation with respect to the other runs. This is due to 
the improved resolution of these simulations that can resolve higher 
densities at earlier times. The star formation rate of the SWio,448 
at lower redshift agrees very well (in fact the black and the grey 
curves are nearly indistinguishable) with that of the SW simula- 
tion, which has the same box size but less particles. This confirms 
the numerical convergence of our simulations. 

The behaviour shown in the star formation rate is very similar 
in nature to that of the total neutral hydrogen evolution of Figure 
[3] (right panel). Here we plot Q.m, which is defined as the contri- 
bution of neutral hydrogen to the total critical density (the values 
have been multiplied by 1000). The neutral hydrogen fraction is 
followed in the simulation self consistently with the assumed aver- 
age UV background and not rescaled a-posteriori using a different 
UV background. Note again that the neutral hydrogen content of 
the WW simulation is about a factor 5 higher than that of the SW 
and this is due to the fact that in the WW simulation the gas is 
colder and more concentrated in the potential wells of galaxies and 
thereby is significantly more neutral. Even in this case the SW2o,5i2 
and SW5,32o simulations embrace all the SW runs as for the star 
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Figure 4. Evolution of the total Hciy (upper panels) and the total f^ovi (lower panels) as a function of redshift for some of the hydrodynamical simulations 
of Table I. Left Panels: effect of different wind strengths and implementations. Middle Panels: effect of different IMFs. Right Panels: resolution tests. 



formation rate plot. The SWkr and SWay agree well with the SW 
confirming that feedback and resolution/box-size effects have more 
impact on the simulations than effects due to the particular choice 
of the IMF as long as the stellar spectra are not self consistently 
taken into account. The difference between WW and SW simula- 
tions is somewhat more pronounced than that found by Nagamine, 
|et al.|p004^ and this is probably due to the fact that our simula- 
tions include metal cooling, thereby increasing the neutral hydro- 
gen amount in the haloes (we will come back to this point later). 
The MDW is in between SW and WW and again the SWio,448 is in 
good agreement with the SW, demonstrating that, at least for these 
quantities, numerical convergence has been reached. 

At the end of Section [43] we compute the DLAs contribution 
to the total neutral hydrogen content (JIdla) and we compare it 
with the results of|Pontzen et ah] ([2008} and the observational esti- 
mates of jProchaska et al.|^2005) and jPeroux et al.| ( [2005| l. 

As a further check, in Figure [4] we plot the evolution of the 
contribution to the total density of two of the ions that trace the 
high redshift IGM at relatively low density and for which there are 
some observational constraints l |Schaye et al.|2003{ [Aguirre et al.| 
2007|: CIV (AA 1548.204, 1550.781 A) in the upper panels and 



OVI (AA 1031.927, 1037.616 A) in the lower panels. Although in 
the rest of the paper we will focus on the distribution of metals 
around DLA systems, we would like to address briefly the evolu- 
tion of ionization species in the IGM as a whole. We use the redshift 
outputs at 2: = 2.25, 3, 4 and we extract the values for ficiv and 
Qovi summing over all the gas particles. In order to obtain the ion- 
ization fraction for the two elements it is necessary to multiply the 
abundance of a given metal, carried by each particle, by its ioniza- 
tion fraction that depends on density and temperature. We use the 
CLOUDY code (Ferla nd et aL|1998} to compute a-posteriori the rel- 
evant fractions for each gas particles. We choose the HMOS option 
in CLOUDY, which consists of a UVB made by QSOs and galaxies 
with a 10% photon escape fraction and which is in agreement with 
other observational constraints ( [Bolton et al.|2005y 

In Figure[4]we explicitly show how different effects impact on 
Qciv and Q,ovi '■ in the left panels we compare simulation with 
different wind strength and implementation (SW, WW and MDW); 
in the middle panels we test the effect of different IMF (SW, SWkr 
and SWay); in the right panels we present resolution tests for all 
the SW runs (SW, SW5,320, SWio,256, SWio,448 and SW2o,5i2). 
Here we discuss each panel of the figure: 
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• Left Panels (feedback physics): for the SW simulations the to- 
tal amount of the two ions increase by a factor of about 3 between 
z = 4 and z — 2.25. The WW simulation contain ~ 2.5 times 
more OVI and CIV compared to the SW case. The reason is that 
with such weak winds the metals remain very close to the dense en- 
vironments around galaxies and cannot reach the low-density IGM 
in contrast with the SW case (see for example [Theuns et al.|2002| l. 
In these regions close to galaxies the ionization fractions are tipi- 
cally larger than in the voids. The MDW CIV and OVI evolution 
is similar to the energy driven implementations although it has a 
higher normalization than SW and WW: this is due to the fact that 
in the former implementation the winds are more efficient in enrich- 
ing the IGM around the mean density with metals and this produces 
larger values of Slciv and Qovi reached at z = 2.25. 

• Middle Panels (IMF): The trend in the redshift dependence of 
f2civ and fiovi is the same for SW, SWkr and SWay, while there 
are differences in the normalization. The Kroupa and Arimoto- 
Yoshii IMFs result in values for the CIV and OVI density that are 
respectively -^1.5 and ~3 times higher than the standard Salpeter 
case. This is due to the fact that Salpeter IMF results in an excess 
of low mass stars, Kroupa IMF produces a smaller number of mas- 
sive stars than the other two but twice as many stars in the range 
of mass 0.3Mq < m < 3AIq, while with the Arimoto-Yoshii (or 
"top-heavy") IMF there is a larger contribution from massive stars 
(see Figure 3 in |Tomatore et al.|2007| ). For these reasons, the SWay 
run produces more oxygen and carbon than the other simulations, 
while the SW run (with Salpeter IMF) is less efficient in producing 
these metal species and SWkr is in between the other two. Inter- 
estingly the MDW run with Salpeter IMF (triple-dot-dashed cyan 
lines in the left panels) produces more CIV than the SWay, while 
the amount of OVI is comparable, confirming the high efficiency of 
momentum driven winds in enriching the IGM. 

• Right Panels (resolution tests): the aim of these two panels is 
to show that our analysis is robust against box-size and resolution 
effects. We plot Jlciv and f2ovi for all the SW runs and it is clearly 
evident the convergence of the results (especially for the Qovi and 
the SWio,448 simulation). As we already found for the SFR and 
the ilm, the SW2o,5i2 and SW5,32o simulations embrace all the 
SW runs. 

A comparison of these findings with those by |Oppenheimer| 
|& Pavel ( |2006| [2008 1 shows that the evolution of ficiv is some- 
what faster and the normalization higher than what they find. In 
particular our MDW is even more discrepant than SW and WW 
with respect to the analogous momentum-driven winds simulation 
("mzw") by jOppenheimer & Dave| ( |2006| l. However, it is difficult 
to compare properly given the different resolutions, box-sizes, 
feedback implementations and details on chemical evolution model 
used. 

In the last part of this Section we briefly discuss the other sim- 
ulations described in Table I. The simulations SW, SWio,448 and 
SWio,256 are characterised by the same box size but different num- 
ber of particles, so that they allow us to carry out a test of stability of 
our results against numerical resolution. As already mentioned, for 
the SWio,448 the star formation rate trend at high redshift is simi- 
lar to that of SW5,32o (the simulation with the highest resolution), 
while it agrees at lower redshift with SW, which has the same box 
size. Instead the SWio,256 has lower SFR (not plotted) than SW at 
high redshift because, for a given box-size, a lower number of par- 
ticles means a lower number of structure resolved in the simulation 
and thereby a lower SFR. As for the evolution of JIhi, f^civ and 




10' 10' 10' 10^" 10" 10^^ 
Mass [h"' Mq] 



Figure 5. Differential halo mass function (number of haloes of a given mass 
per unit mass [h~^ Mq] and per unit volume [Mpc^/h^]) for runs SW, 
SW5.320, SW20,512 at redshift z = 3, compared with |Sheth & Tormen| 
I .1999j prediction. 

f2ovi, the SWio,448 follows the results of the SW almost perfectly. 
This is also true for other analysis made in this paper (phase dia- 
grams, haloes properties, etc.), thus further confirming the numeri- 
cal convergence of our results on these simulations properties. The 
run SWwDM has a different initial power spectrum P{k) respect all 
others simulations: power is suppressed at small scales due to free 
streaming of the dark matter particles that in this case are "warm" 
and have mass equal to 1.2 keV (corresponding to a suppression 
scale of around 50 comoving h^^kpc). As a consequence, the star 
formation rate density at high redshift is significantly reduced with 
respect to the CDM case, with a correspondingly smaller number of 
haloes having mass comparable or smaller than the free-streaming 
mass scale. We run this simulation in order to test whether the mod- 
ified P{k) could produce haloes that better fit the column density 
distribution function for small values of A^hi, that will be discussed 
in Sectionl431 



4 THE NEUTRAL HYDROGEN DISTRIBUTION 

In this Section we focus on the properties of the neutral hydrogen 
around the simulated galactic haloes. The aim is to study DLAs 
properties and statistics extracted from our simulations and to com- 
pare them with the latest observational data and the simulation re- 
sults from other groups l |Nagamine et al.|2004[|2007{|Pontzen et al.| 
[2OO8I 1. 

4.1 Identifying haloes 

We identify haloes in the simulations by running a parallel Friends- 
of-Friends (FoF) algorithm with a linking length which is 0.2 times 
the dark matter mean interparticle spacing. The number of haloes 
of a given mass per unit mass and per unit volume as a function of 
their mass is shown in Figure|5]for the SW, SW5,32o and SW2o,5i2 
runs at redshift z = 3, along with the |Sheth & Tormen| ( |I999[ > mass 
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function, so as to test at the same time thie effect of resolution and 
of box size. Again, the SW2o,5i2 and SW5,32o embrace the refer- 
ence case SW: the mass function of SW5,32o extends to small halo 
masses due to its better resolution, while SW2o,5i2 produce more 
massive haloes due to its larger box size. In our largest simulation 
we have three haloes above W^^h~^ Mq and about 1000 haloes of 
masses above W^'^h~^ Mq at z = 3. 

We follow the analysis made by |Nagamine et al7( 20O4, 2007^ 
to realize a mock DLA sample: after having identified the haloes 
and their center of mass, we interpolate with a TSC (Triangular 
Shaped Cloud) algorithm the comoving neutral hydrogen mass den- 
sity around the center of mass of each halo on a cubic grid; then we 
'collapse' the grid along a random direction and we obtain a set of 
neutral hydrogen column densities for each halo. Thus the column 
density reads: 

AThi = SiPi,Hie/mp(l + z)' , (5) 

with nip the proton mass and e = i/wgrid the linear dimension of 
the single grid cell. Here / is the size of the box around the halo and 
rigrid is the number of grid points. 

Tipically, for the most massive haloes, we use cubes of size 
200 comoving /i^^kpc with 32'^ grid points (e = 6.25/i^^kpc). In 
such a way, we increase the DLA total redshift path and we sample 
32^ X A^haiocs HI column densities along lines-of-sight per simu- 
lated box. 

We have carried out some tests changing the number of grid 
points in order to study the effect of the sampling size on the neu- 
tral hydrogen distribution. The choice of the grid points number is 
crucial because too many points produce a sampling size below the 
resolution of the simulation and consequently an "oversampling" 
of the HI mass density with large statistical fluctuations, while too 
few points produce a smooth statistic which could not be repre- 
sentative of the real density field. At the end we found that the 
best compromise was to use 32'^ grid points, corresponding to e 
^ 4.5 X softening length (which is also the typical value of the 
SPH smoothing length in the outskirts of the haloes), differently 
from |Nagamine et al.| ( |2004) , who choose instead e ~ softening 
length. 

Figure |6] shows the HI column density maps extracted as ex- 
plained above, but with a finer grid subdivision of 256 points, for 
the same massive halo in the WW (upper panel), SW (middle panel) 
and MDW (lower panel) runs at z = 3. The HI density at each 
pixel has been projected along the line-of-sight in the z direction. 
In the figure it is visible the effect of the winds: in the WW run 
(upper panel) high column density gas is more concentrated inside 
the central halo and inside some substructures. Also the column 
density values reached are higher than for the other runs. In the SW 
run (middle panel), the gas is more spread around the central haloes 
and the substructures. Finally if we consider HI column densities 
above the DLA limit of A^hi = 2x10^° cm^^ we see that the cen- 
tral halo in the MDW run has the largest cross-section. We discuss 
further about this point in Section [43] 

4.2 Properties of the haloes 

For the haloes identified with the FoF algorithm we compute mean 
quantities that could be relevant for the following analysis, such as 
the star formation rate inside each halo, the mass-weighted mean 
total metallicity and the mass-weighted mean neutral fraction of 
hydrogen (HI/H). We plot in Figures |7][8]and|9] our findings at red- 
shift z — 3, only for haloes having mass greater than 2 x 10*/i~^ 
Mq, that are resolved with at least 100 dark matter particles. The 
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Figure 6. HI column density maps in a slice around the same massive halo 
in the WW (upper panel), SW (middle panel) and MDW (lower panel) runs 
at 2 = 3. 
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Figure 7. Star fomiation rates (in Mq per year) plotted as a function of halo mass for SW, WW and MDW runs at 2 = 3. 
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Figure 8. Mean total metallicity in solar unit plotted as a function of halo mass for SW, WW and MDW runs at 2 = 3. 



properties of the haloes below l(fh~^ Mq, resolved with ~ 1000 
particles, should not be trusted at a quantitative level, but we prefer 
to show them in order to appreciate the increase in the scatter at 
low masses. 

For the most massive haloes, the WW simulation shows very 
high star formation rates that are even a factor 100 larger than that 
of the corresponding haloes in the SW run: this is expected due to 
the smaller efficiency of the weak winds in suppressing star forma- 



tion (see the left panel of Figure [3](. The MDW run reaches nearly 
the same SFR values of the SW for the most massive haloes, while 
for intermediate mass haloes (10^'^~"/i"^ A/q) the MDW haloes 
have higher SFR than SW ones. For the smallest haloes of around 
IQ^h"^ Mq the SW run shows a slightly larger scatter in the SFR 
values. Overall, the bulk of the haloes, of masses between 10^ and 
lO^^/i"^ Mq, that are likely to host DLA systems, have star for- 
mation rates of about 0.1 M0/year. For haloes above 10^*^ Mq 
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Figure 9. Mean neutral fraction of hydrogen (HI/H) plotted as a function of halo mass for SW, WW and MDW runs at 2: = 3. 



the star formation rates are different between the models and usu- 
ally around 1-30 MQ/year. These star formation rates are in general 
agreement with those of the population of faint Lyman-a emitters 
found recently by "Ranch et al.| ( [2b08 1, in which a link between the 
Lyman-Q emitters and the DLAs is suggested. 

For the same reasons metallicities of haloes in the WW run 
are higher than in the SW and MDW simulations (see Figure [Sjl. 
Comparing this result with Figure |4] one can derive the following 
general picture: in the WW run metals remain locked inside the 
haloes (or stay close to the haloes), while in SW and MDW they 
are able to reach the IGM and enrich it: this is particularly true 
for the less massive haloes in the MDW. In fact, the MDW small 
mass haloes have in general metallicity lower than the correspond- 
ing SW and WW haloes. The wind implementation is thereby very 
effective in devoiding the galaxies of star forming cold gas that is 
enriched at a level of about 0. 1 Zq . The most important result of 
this panel is however that the different wind implementations show 
different mass-metallicity relations especially for haloes of masses 
above lO^'^ft~^M0: while WW and MDW show a correlation (al- 
though with different amplitude), the SW results seem to produce 
very little correlation between mass and metallicity or possibly a 
weak anti-correlation. This different metallicity pattern in haloes 
could be important when compared with observation (e.g. |Maiolino| 
|et al.|2008 ) and could possibly allow to discriminate between dif- 
ferent feedback scenarios. 

Finally we note that the largest metallicities appear to occur in 
haloes of small masses. There are two different motivations for this 
to happen. The first motivation is numerical: as we mentioned at 
the beginning of this Section, both these small haloes and the wind 
model are not well resolved, thereby their large metallicities are, at 
least partially, the result of some numerical artifact. Furthermore, 
we also checked that there is a physical motivation: these small 
haloes are strongly affected both by self enrichment at early epochs 
and by the metal enriched winds blowing from bigger haloes close 



Redshift 


Run 


slope a 


P 




SW 


0.57 


3.83 


z = 4 


WW 


0.46 


3.79 




MDW 


0.69 


4.11 




SW 


0.77 


3.75 


z = 3 


WW 


0.49 


3.61 




MDW 


0.85 


4.06 




SW 


0.62 


2.99 


2 = 2.25 


WW 


0.52 


3.45 




MDW 


0.92 


3.88 



Table 2. Fitting parameters a and /3 of Eq. j6) for runs SW. WW and MDW 
at redshit z = 4, 3 and 2.25. 



to them. Addressing quantitatively these issues would require more 
numerical work, which is beyond the scope of this work. 

In Figure|9]we show the mass-weighted mean neutral fraction 
of hydrogen (Hl/H) inside haloes, a quantity which is closely re- 
lated to the DLA properties. The general trend reflects those of the 
right panel of Figure |3] for a given halo mass the WW simulation 
displays the highest value, while the SW and the MDW have com- 
parable values. In all the panels, there are present two different sets 
of values: above and below HI/H ~0.01. The first is associated to 
the most massive haloes and the second to the least massive ones. 
This is due to the fact that the most massive haloes contains many 
particles above the density threshold pth for which the neutral hy- 
drogen fraction is equal to the fraction of mass in cold clouds /c 
(see Eqs. (1) and (2) in Section 2), boosting in this way the neu- 
tral hydrogen content of these haloes. Less massive haloes have 
instead neutral fractions that are set by the physical conditions of 
the gas and by the UVB background and are typically more ionized 
(/hi /c) than the most massive haloes. 
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4.3 The DLA cross-section 

For each halo of total mass Mtot we derive the DLA cross-sections, 
o"DLA (in comoving units), by selecting and summing up the area 
of all the cells with a column density (determined as explained in 



of the form: 



log (TOLA = a(log A/tot - 12) + /? . 



(6) 



for all the simulations made. Results are shown in Figure [TO] 
where we plot the cross-sections as extracted from the SW at 2; = 
2.25, 3, 4 with overplotted fits for SW itself (dashed black line), 
WW (blue dot-dashed line) and MDW (cyan triple-dot-dashed 
line). Fitting parameters a and 13 for each of these runs are show 
in Table II. The continuous red lines in the middle panel fit the up- 
per and lower envelopes of the SW distribution at 2 = 3 and were 
drawn with the aim of highlighting the scatter in the distribution 
of cross-sections expecially at low masses. The parameters of the 
upper envelope line are a — 0.5 and /3 = 3.6 (very similar to the 
WW case), while those of the lower envelope line are a = 1.2 and 
P = 4.02. In the next section we discuss the effect of taking into 
account this scatter in a very conservative way. 

It is clear that for the WW simulation the cross-sections are 
larger than in the SW case for the small mass haloes: this is due 
to the fact that the wind is more effective in expelling cold gas 
from small mass haloes than from the most massive ones. For the 
most massive haloes the content of cold gas is much larger and 
this trend is inverted: the WW run has a smaller cross-section than 
in the SW and MDW cases because the cold gas is located very 
close to the dense galactic environments. The MDW trend reflects 
the fact that momentum-driven winds are more directly related to 
the halo properties: smaller mass haloes have greater loading factor 
and winds become very efficient as in SW case. On the contrary, 
for massive haloes the MD winds become less efficient: the cold 
blob of gas around the halo has a more regular shape than in the 
SW case and determines a larger cross-section than the SW. To 
sum up, as one can see from Figure 6, the qualitative trends are 
the following: for WW the gas is concentrated in the center of the 
haloes; the SW expels gas in a quite violent way and the shells 
of expanding material fragment; the MDW blows less gas but in a 
more homogeneous way than SW. 

A comparison of the results from the SW run with those ob- 
tained by |Nagamine et alT]j2004[[2007| l at similar resolution shows 
that we obtain on average a 20 per cent shallower slopes and 5% 
smaller (3 parameters: this is likely to be due to the faster wind 
speed adopted for the strong wind case (600 km/s vs. 484 km/s) that 
will probably slightly reduce the normalization (3 devoiding galax- 
ies of cold gas. Furthermore, the fact that our simulations have the 
metal cooling implemented could result in a larger cross-section 
for smaller haloes because the quantity of cold gas increases, while 
for the most massive haloes the metal cooling enhances the amount 
of stars at the expenses of cold gas. These two effects determine 
a shallower slope for the cross-section fitting. In the WW case the 
run should be compared with the P3 run of |Nagamine et al.]j2004| 
|2007| even though the wind speed adopted in this study is 100 km/s 
vs. 242 km/s and the mass resolution is about 10 times better, but 
even in this case the trend is confirmed and we find a shallower 
slope, of about 40 per cent, and a smaller normalization value (by 
~ 10 per cent). 

In Figure [To] it is clearly visible, especially at redshift z — 
2.25, a discretization in the cross-section values, more pronounced 
than in |Nagamine et al.| ( [2004[|2007[ l. This is due to our final choice 



of the linear dimension of the single grid cell e which is somewhat 
larger than that adopted by |Nagamine et al^j2004^ . 



4.4 The incidence rate of DLA systems 

Having obtained the mean relation for the DLA cross-section as a 
function of halo mass from the previous Section, we are now able 
to calculate the cumulative number of DLAs per unit redshift (or 
rate of incidence) using the equation: 



dz 



'-{> M, 



dr 
dz 



ni,{M', z) (Jdla(M', z) dM' , 



(7) 

where nh{M, z) is the Sheth & Tormen 1 1999^ dark ma tter halo 
mass function and dr/dz = c/i/o\/f2m(l + 2)'' + i^A. Follow- 
ing Nagamine et al. (2004), we use this equation in order not to 
be sensitive to dark matter haloes with masses below the resolu- 
tion limit of the simulation. This is a common problem when one 
tries to compute the number density of DLAs based on a cosmo- 
logical simulation that does not resolve all small mass haloes that 
may host a DLA: the underlying assumption in the rest of the pa- 
per is that haloes below ICp Mq are not able to produce DLA 
systems. Moreover small box size simulations cannot produce very 
massive haloes. To overcome these limitations in Eq. (jTJ we con- 
volved the Shet h & Tormen| ( |1999[ l halo mass function with the 
measured relationship between DLA cross-section and halo mass, 
thereby correcting for incompleteness in the resolved halo abun- 
dance of our simulations. In doing this we extrapolate the power- 
law fit of Eq. l|6j both at high and low masses. This relation presents 
an intrinsic scatter: looking at Figure [lO] one can note that a given 
o"DLA value corresponds to different halo masses, especially at low 
cross-sections. As mentioned in the previous section we check the 
impact of this scatter fitting the upper and lower envelopes of the 
cross-sections distribution of the SW run at z = 3. 

In Figure [11] we show dNoLA/dz for SW, WW and MDW 
runs at z = 3. The three different runs have identical initial power 
spectrum and cosmological parameters so they all have the same 
theoretical dark matter halo mass function. Therefore the differ- 
ences in Figure[TT|reflects what we found in previuos Section about 
DLAs cross-sections (see in particular the central panel in Figure 
|10^ . In fact at low masses WW simulation produce haloes with 
cross-sections ~3 times higher than SW and MDW, and corre- 
spondingly dNuLA/dz reaches higher values. At greater masses 
the trend is inverted and WW curve stays below the other two. 
As we expected dNohA/dz curve for MDW simulation is al- 
ways above the SW one and well above WW one in the high 
mass tail. The shaded region shows the observational estimate 
log(dAi'/dz)=— 0.60±0. 10 at z = 3 recovered by |Nagamine et af] 
l |2007) from the obs ervational results based on SDSS QSO spectra 
of Prochaska et al.|f2005| l. Finally, the continuous red lines show 
the results considering the scatter in the cross-section vs mass rela- 
tion for the SW run. In one case (where the curve is almost equal 
to that of the WW case) one would need only haloes with mass 
greater than ~ IO^ '^/i'^Mq to fit the dNuhA/dz statistic, while 
in the other case one should go down to much less massive haloes, 
that are not well resolved by our simulations. 

Our findings for the DLA abundance per unit redshift are 
slightly different than those obtained by Nagamine et al.| ( |2007| l: 
moving to low masses their distribution tend to flatten, while ours 
are somewhat steeper and this is due both to the different val- 
ues of the cross-section, the different dimension of the grid cells 
used, and the different cosmological parameters for the linear dark 
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Figure 10. Triangles represent the DLA cross-sections (in comoving units) as a function of the total halo mass for the SW run at 2 = 2.25,3,4. The overplotted 
lines show power-law fits, of the form log ctjjla = "^(log A/tot — 12) + /3, to the data points for the SW itself (dashed black line), WW (dot-dashed blue 
line) and MDW (triple-dot-dashed cyan line) runs. The continuous red lines in the middle panel fit the upper and lower envelopes of the SW distribution at 

2 = 3. 



matter power spectrum (amplitude and slope). As a result, at low 
masses, the dNoLA/dz values in our simulations are about a fac- 
tor ~ 2 greater than the corresponding ones of [Nagamine et aE] 



( 2007 I, while at high mass the trend is inverted and our dNuhA/dz 
decreases somewhat more rapidly than Nagamine et al. results. 
To summarize, our results reproduce the observational data of 
Prochaska et al.| ( [2005 [ l slightly better than those of [Nagamiriel 



et al. 



1 2007 1 which slightly underpredict dNuLA/dz in the range 



of masses considered. 

As far as the redshift evolution of DLA rate of incidence is 
concerned, the behaviour of the three different simulations reflects 
that of Figure [lO] For the sake of brevity we do not show all the 
plots, but at 2 = 4 the decreasing in dNoLA/dz at high masses 
is much more pronounced than at 2: = 3, due to a lack in massive 
haloes that have not already formed. At 2 = 2.25 both WW and 
MDW curves stay always above the SW one and are more flattened 
than before down to low masses. The overall trend with redshift is 
the same of the left panel in Figure 6 of |Nagamine et al.] ( [2004| >. 

4.5 The column density distribution function 

In this Section we investigate the column density distribution 
function for DLAs usually plotted in the form of f{N), where 
f{N, X)dNdX is the number of DLAs with HI column density 
in the range [A'^, N + dN] and absorption distances in the inter- 
val [X, X + dX\. The absorption distance is given by X{z) — 
S;{l + z'fH,/H{z')dz'. 



In the left panel of Figure 12 we show f{N) for the SW 



(dashed black line), WW (dot-dashed blue line), SW2o,5i2 (dotted 
red line) and MDW (triple-dot-dashed cyan line) runs at 2 = 3. The 
data from SDSS are overplotted along with the F-function fit (black 
diamonds and solid line) of Procha ska et ah] ^2005^ . SW2o,5i2 is 
plotted to test the effect of different resolution and box size and its 
trend is similar to SW. Both are in quite good agreement with data 
especially at low and high column densities. Instead WW overpre- 
dicts the distribution function at large A^'hi as was already found 
by |Nagamine et ar] ( |2004[ l. Finally, MDW agrees extremely well 




Figure 11. Cumulative abundance of DLAs per unit redshift as a function of 
total halo mass for the SW (dashed black line), WW (dot-dashed blue line) 
and MDW (triple-dot-dashed cyan line) runs at redshift 2 = 3. The continu- 
ous red lines show the results considering the scatter (in a conservative way) 
in the cross-section vs mass relation for the SW run. The red shaded re- 
gion indicates the observed cumulative DLA abundance of |Prochaska et al.] 
^2005|l from SDSS data. 



with the data down to the smallest A^hi values. Comparing with the 
results of Naga mine"etaL| pOOil [2007] ) we find a better agreement 
with the observational data both for the SW and MDW, especially 
at low column density and at 2 = 3. 

In the right panel Figure[T2]we split the contribution of the dif- 
ferent haloes to the total f{N) for the SW2o,5i2 at 2 = 3. We use 
the SW2o,5i2 run because the largest box size improves the statis- 
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Figure 12. Left Panel: HI column density distribution function at 2 = 3 for tlie SW (daslied black line), WW (dot-dashed blue line), SW2o,5i2 (dotted red 
line) and MDW (triple-dot-dashed cyan line) runs. Right Panel: Contribution to the HI column density distribution function from haloes of different mass in 
the SW2o,5i2 run at 2: = 3. Cyan triple dot-dashed line refers to haloes with mass in the range 10^ — \0^^h~^ Mq, orange dot-dashed line refers to haloes 
with mass in the range 10^'^ — 10^^h~^ Mq, red dotted line refers to haloes with mass greater than \0^^h~^ Mq while black dashed line refers to all haloes. 
In both the panels the overplotted black diamonds and the solid line show the data points and the fit of |Prochaska et al.]j2005) . 




Figure 13. Left Panel: HI column density distribution function at z = 2.25 for the SW (dashed black line), WW (dot-dashed blue line) and MDW (triple-dot- 
dashed cyan line) runs. Overplotted black diamonds and solid line show the data points and the fit of [Prochaska et al. (2005 1. Run SW20.512 is not present 
here because this simulation was stopped at 2: = 3. Right Panel: the same as the left panel but at redshift z = 4. Dotted red line refers to the SW2o,5i2 run. 



tics on the haloes. Looking at the figure one can easily see that 
haloes with masses above 10^^ h^^ A4q (dotted red line) contribute 
primarily to the large column densities, while the smallest haloes 
with masses in the range 10^-10^'^ h^^ Mq (triple-dot-dashed cyan 
line) contribute significantly to the lower column densities, below 
A^'hi ~ 10^"'* cm"'^. This is not surprising since on average more 
massive haloes have higher gas densities and corrispondingly pro- 
duce larger HI column densities values. 



In Figure [T3] we show the same as Figure [12] but at redshift 
z — 2.25 (left panel) and z — 4 (right panel), to check the redshift 
evolution of the column density distribution function. In the left 
panel, the results for SW20.512 are not shown because this simula- 
tion ended at 2 = 3. SW fits well the observational data for HI col- 
umn densities larger than 10^^ cm~^, but at lower column densities 
there is a discrepancy of about a factor 3. WW and MDW overpro- 
duce f{N) at high A'^hi, but there is a better agreement with data 
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at logiVHi < 20.7, especially for the MDW run. At redshift 2 = 4 
(right panel) SW, SW2o,5i2 and also MDW match very well the 
data while WW does not. The fact that MDW fits the data as well 
as SW at 2 = 4 is due to the fact that high redshift momentum- 
driven winds behave quite similarly to energy-driven ones, while at 
lower redshift the two models tend to differ significantly in terms of 
velocities and loading factors. Thus, we conclude that the different 
wind implementations of the galactic wind feedback show distinct 
predictions for the redshift evolution, and in general the differences 
become larger when moving to 2 ~ 2. At the end the feedback 
model that reproduce better the data is the MDW model. 

At the end of Section [3^ we briefly discussed the SWwdm 
simulation. We run this simulation motivated by the poorer fit to 
observational data for the systems of column densities between 
-|^q20.3-20.8 ^^-2 jjj^j jjjjjjjj in |Nagamine et al.||2004||2007| 

and |Pontzen et al.|^2008[ l. Our reference runs (especially the SW) 
show a similar flattening trend and underproduce the number of 
these systems, by a smaller amount when compared to Naga mine| 
|et al.|p007^ . The idea is that using a different linear power sprec- 
trum we could possibly modify the halo mass function and this 
could impact on the column density distribution function as well. 
We perform the same analysis of the other runs for SWwdm, but 
we find no statistical significant difference for the whole column 
density range. The smaller number of dark matter haloes, com- 
pared to ACDM cosmology, is thus compensated by an increase 
in the cross-section, since in WDM the haloes are usually less con- 
centrated. Thereby we decided not to consider this run anymore 
and we conclude that it appears unlikely that this statistic could be 
better fit by invoking modifications of the linear dark matter power 
spectrum such as warm dark matter. 

Finally, we compute the total neutral gas mass in DLAs using: 

nDLA(2)= ° / f[Nm,X)NmdNm, (8) 

c/hi Pc,0 J1020.3 

where mp is the proton mass, /hi is the neutral hydrogen fraction 
of the gas and pc,o is the critical density at redshift 2 = 0. The 
integration limit goes from 10^" "^ cm~^ (the lower column density 
limit for a system to be identified as DLA) to TVmax = 10^ 
cm~^, this latter chosen to compare with the results of Pontzen 
[etaLlpOOSj . 

In Figure [T4]we show the evolution with redshift of Odla for 
the SW (black solid line), WW (blue dotted fine) and MDW (cyan 



triple-dot-dashed line) runs, along with the result of [Pontzen et al. 



(20081 (green cross) and the observational estimate of Prochaska 



[et al., ( ,2005) (red diamonds) and |Peroux et al.||2005) (orange trian- 
gles). The different simulations' trends reflect those plotted in the 
right panel of Figure [3] the WW run produces the largest SIdla 
value while SW the smallest. The amount of neutral hydrogen in 
DLAs is about a factor two smaller than the total neutral hydro- 
gen in the simulate d volume. The MDW i s in good agreement with 
the value found by [Pontzen et al.| ([2OO8) (1.0 x 10"^ in the red- 



shift range 2.8 < 2 < 3.5) and with the data by |Peroux et al.| 
( |2005 ), while at redshift below 2 ~ 3 it is slightly larger than the 
recent measurement made by [Prochaska et ar] ( [2005[ > using DLAs 
in SDSS spectra. 

So far we did not comment much about self-shielding effects. 
A full treatment of the self-shielding would require radiative trans- 
fer and this again would be an approximate (usually a-posteriori) 
scheme. We decided to rely on the multiphase ISM model to im- 
plictly account for self-shielding motivated by the two following 
facts: i) many observational properties of DLAs have been repro- 
duced in such a way ( |Katz et al.|1996i|Nagamine et al.|2004|[2007l l; 
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Figure 14. Redshift evolution of Hdla for the SW (black solid line), WW 
(blue dotted line) and MDW (cyan triple-dot-dashed line) runs. Overplotted 
are the result of Pontzen et al. (20081 (green cross) and the observational 
estimate oPProcha ska et al. 1 2005) (red diamonds) and ,Peroux et al.,^2005) 
(orange triangles). 



ii) the recent results by Pontzen et al. 1 2008 1 seem to suggest that 



even when applying crude radiative transfer approximations the 
main properties of DLAs do not change significantly. However, to 
better check the impact of this criterion, we also apply a further ap- 
proximation that is useful in order to decouple the star formation 
from the multiphase prescription and assume that the gas particles 
above densities of nn = 10~^ cm""' are fully neutral (following 
[Haehnelt et al.|1998) . We recomputed the column density distribu- 
tion function for a few simulations and found that the differences 
are not large. All the column density distribution function increase 
by a fixed overall amount of about 0.2 dex and are still in broad 
agreement with the observations. Since we believe that our multi- 
phase ISM model is more refined and more physically motivated 
than this criterion we decided not to present results for this second 
simpler assumption. 



5 SIMULATED QSO SPECTRA 

For each simulation performed we have extracted several physi- 
cal quantities interpolated along lines-of-sight (LOSs) through the 
box, following the procedure of |Theuns et al.|jl998| ) (Appendix 
A4). The optical depths of the simulated QSO spectra are drawn 
in redshift space taking into account the effect of the IGM peculiar 
velocities along the line-of-sight, Wpec.n- As already mentioned in 
Section 2, our simulations follow self consistently the evolution of 
H, He, C, O, Mg, S, Si and Fe. Using the CLOUDY code < Ferland| 
|et al.|[l"998] l, we then determine the ionization fractions for some 
ions that trace the high redshift IGM: CIV (AA 1548.204, 1550.781 
A), OVI (AA 1031.927, 1037.616 A) and Sill (A 1526.707 A), with 
the possibility to extend the analysis to many others. The simulated 
flux of a given ion transition at the redshift-space coordinate u (in 
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Figure 15. Comparison of physical quantities extracted from tliree simulated LOSs passing through the center of mass of the second most massive halo for the 
SW (solid black line), WW (dashed blue line) and MDW (dotted red line) runs at z = 3. From top to bottom the figure shows: total gas density [nn/cm'^], HI 
density [njj/cm^], total metallicity (in solar unit), temperature [K] and peculiar velocity field [km/s]. 



km/s) is F{u) — exp[— r(u)] with: 

t{u) — '^'^ / dxnx{x) V \u — X — v^^^\^{x), b{x) \ dx , 

(9) 

where (jq.x is the cross-section of the particular ion transition, 
H{z) is the Hubble constant at redshift z, x is the real-space coordi- 
nate (in km s~^), b = {2kBT / my^c?)^ ^"^ is the velocity dispersion 
in units of c, V is the Voigt profile. 

In Figure[T5]we compare three LOSs along the center of mass 



(CM) of the second most massive halo in the SW, WW and MDW 
simulations, at redshift 2 = 3. We interpolate physical quantities 
along the LOS using the SPH kernel of each particle. The SW is 
represented by the continuous black line, the WW is represented 
by the dashed blue line and the MDW by the dotted red line. From 
top to bottom we show the gas density, neutral hydrogen density, 
total metallicity, temperature and peculiar velocity. The a;— axis is 
in km/s and represents at a given redshift the size of the box in 
real space. Looking at the first two top panels it is visible a peak 
both in the gas and in the neutral hydrogen densities: this corre- 
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Figure 16. Left Panel: DLAs velocity width distribution in runs SW20,512 (dotted red line), WW (dot-dashed blue line), MDW (triple-dot-dashed cyan line) 
and SW2o,5i2 smoothed (solid black line) at z = 3. Right Panel: DLAs velocity width distribution in runs SW2o,5i2-''"Hl (solid red line), WW_thi (solid 
blue line), MDW.tjji (solid cyan line) and SW2o,5i2 smoothed (solid black line) at z = 3. In both the panels the overplotted black diamonds show the 
observational data of |Prochaska et al!]j2008) . 



spends to the CM of the halo. All simulations have more or less 
the same structure for the peak in density but in the MDW case 
(red dotted line) it is a bit more extended than the other two. The 
metallicity panel clearly shows the effect of the different winds im- 
plementations and velocities. The SW run produces a distribution 
of metallicity more extended than MDW and WW while the WW 
distribution is narrower than the MDW one. This confirms the fact 
that for massive haloes MDW is less efficient than SW in extract- 
ing metals and distributing them in the surrounding IGM (because 
of the high velocity dispersion of massive haloes that results in a 
small loading factor for the wind particles even if the correspond- 
ing wind velocities are high). The temperature panel presents the 
same trends as the metallicity panel: the high temperature peak is 
larger for SW with respect to MDW and WW, demonstrating that 
the hot SW gas can be spread out to large distances. Quite interest- 
ing in this panel is the drop in temperature at the center of the halo 
where the gas density is high and the cooling efficient. The bottom 
panel shows the peculiar velocity field. One can clearly see the ef- 
fect of the expanding wind: a strong discontinuity in the peculiar 
velocity produced by the shock, with a negative peak (shell coming 
toward the observer) and a positive one (shell expanding along the 
LOS in the opposite direction). Again the SW have the larger distri- 
bution and WW the narrower MDW have larger peaks (in absolute 
values) than SW but not as extended, that are produced by a fast, 
but less mass-loaded, wind. We note that the peculiar velocity fields 
account for the observed wind velocities in local starburst galaxies 
and that these velocities do depend on the wind model. 

5.1 Velocity width distribution 

In this Section compare observational data and simulation predic- 
tions on the observed velocity width distribution for low-ionization 
species. We use the quantity Atigo defined as the extent in veloc- 
ity (redshift space) that embraces 90% of the total integrated op- 



tical depth as in jProchaska & Wolfe| ( [1997] l. The procedure to re- 
cover this statistics is the following: for every simulation we take 
one thousands spectra extracted along the CMs of the most mas- 
sive haloes in the box. Then for every spectrum we compute the 
absorption profile of the ion considered (Sill, A 1526.707 A) and 
finally from this we determine the values of Augo associated with 
that particular LOS (i.e with that particular DLA). To deal with the 
fact that Sill occurs in self-shielded regions we use the ISM multi- 
phase model in the same way as we did for the HI. Following Eqs. 
(1) and (2), we assign to each particle in the LOS a mass in Sill: 

msai = /siii msi {p < pth) (10) 
"isiii = /c?nsi {p>Pth), (11) 

where msi is the mass in Si of the particle (determined self con- 
sistently inside the code), fsni is the Sill fraction as determined 
by CLOUDY on the base of the temperature and the density of the 
particle, fc and pth are, respectively, the fraction of mass in cold 
clouds and the star formation threshold (see Section 2). 

It is well known that the observed distribution of velocity 
widths has a median too large to be accomodated within stan- 
dard cold dark matter cosmology jProchaska et al. 2008 ; Pontzen| 
|et al.||2008l l. In the left panel of Figure |16| we show the DLA 
velocity width distributions at redshift z — 3 for our refer- 
ence runs SW2o,5i2 (dotted red line), WW (dot-dashed blue line), 
MDW (triple-dot-dashed cyan line) along with observational data 
of |Prochaska et al7| l |2008^ . The latter consist of a set of 113 mea- 
surements of metallicity, redshift and Avgo taken with HIRES, 
ESI and UVESA'LT spectrographs for DLAs in the redshift range 
1.5 < z < 4.6. We consider this sample as a whole since it does 
not show any redshift evolution |Pontzen et al.|20(j8] l and we com- 
pare with our 2 = 3 outputs. 

It is clear that our reference runs dramatically fail to reproduce 
observational data even to a larger extent than found bj|Pontzen] 
[eTaTlpOOSl l. In fact SW2o,5i2, WW and MDW overpredicts the 
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number of small velocity systems and at the same time they are not 
able to produce systems with velocity width greater than 100 km/s. 
As we expect, MDW does a better job than the other two runs, 
since it is more efficient in expelling metals, but it is yet far away 
from fitting data. Checking one by one a large number of extracted 
spectra we find that for a given LOS there is a significant number 
of small regions in the optical depth, all of them with sizes that are 
too small in redshift space. Basically our simulations, regardless to 
the particular wind implementation, spread around the haloes small 
clumps of enriched materials, that are in a wind phase, and are not 
able to enrich uniformly the surrounding IGM. This of course also 
results in the overproduction of small velocity systems seen in the 
left panel of Figure [T6] 

We checked that even using a simulation with higher resolu- 
tion like SWio,448, which in principle should start enriching the 
IGM earlier, we are not able to reproduce observational data. To 
overcome the problem we follow two other post-processing criteria 
with the goal of finding out some guidelines to improve our future 
simulations: 20. 512 smoothed shown in both the panels of Fig- 
ure [16] as the black solid line and the "thi" series plotted in the 
right panel of Figure [T6] 

For the thi series we take the original SW2o,5i2, WW and 
MDW runs and then we follow [Pontzen et aL| p008) assummg 
that Sill is perfectly coupled to HI so that for solar metallic- 
ity Mx/Mh = 0.0133 and n(SiII)/n(HI) = n(Si)/n(H) = 
3.47 X 10"^ In this way Sill is spread more efficiently around 
and inside the haloes. Moreover here we sum, for a given LOS, all 
the different Awgo contribution due to the small metal clumps try- 
ing to mimic a more uniform enrichment. The results are shown in 
the right panel of Figure [16] SW2o,512-T"hi is the solid red line, 
WW.thi the solid blue line and MDW.thi the solid cyan line. 
SW2o,512-Thi works better than the other two mostly because it has 
a greater number of bigger haloes where the Augo is correspond- 
ingly higher (see below), but as one can see this post-processing 
procedure is not yet enough to match the observations. We find a 
trend similar to what |Pontzen et aL| (2008) showed in their Figure 
9: all the different thi series reproduce the small velocity tail of the 
distribution but fail to produce systems with velocity width larger 
than ~300 km/s. 

Instead, for 20, 512 smoothed we extract simulated spec- 
tra after having increased the smoothing length associated to each 
particle along the LOS and having recomputed the metallicity. Ba- 
sically we set the new smoothing length to 500 h^^ comoving kpc 
and the result is the 20,^12 smoothed curve (black solid). With 
this prescription we can fit very well observational data both at high 
and small velocities. Of course spreading the metals a-posteriori 
over a scale of hundred /i~^kpc is not self consistent but in this 
way we just want to see which phenomenological prescription can 
be used in order to fit the data. 

Either a pre-enrichment of the IGM at higher redshift that 
could pollute the gas particles in a more uniform way (e.g. |Tor-| 
|natore et al.pOOTl l or some missing physical ingredient such as 
smsill scale turbulence that could be effective in mixing metals in a 
more efficient way at galactic scales (e.g. |Scamiapieco & Bruggen] 
|2008[ll^ichino et al.|2008| l, could help in reproducing the observed 
values. From this point of view a smoothing scale of 500 /I'^kpc 
roughly corresponds to a velocity (in redshift space) of ~50 km/s, 
values that are consistent with those presented in |Scannapieco &] 
|Bruggen| ( |2008) . 

Finally, in Figure [TT] we test the contribution to the DLAs ve- 
locity width distribution from haloes of different mass at 2 = 3. 
We use the three models that fit best the Augo statistic, the thi 




100 200 300 400 500 600 700 
Avjd [km/s] 

Figure 17. Contribution to the DLAs velocity width distribution from 
haloes of different mass in runs SW20,512-thi (solid + dashed red lines), 
WW.rjji (soHd + dashed blue lines), MDW.rjjj (solid + dashed cyan lines) 
at 2 = 3. For the different runs, dashed lines refer to haloes with mass lower 
than IO^^'^/i^^Mq and solid lines refers to haloes with mass greater than 
10^^'^h~^ Mq. The overplotted black diamonds show the observational 



data of |Prochaska et al.|j2008) . 



series: SW2o,512-T"hi (solid -1- dashed red lines), WW_rHi (solid + 
dashed blue lines), MDW .thi (solid -1- dashed cyan lines). For the 
different runs, dashed lines refers to haloes with mass lower than 
W^'^'^h~^ Mq and solid lines refers to haloes with mass larger than 
W^'^'^h~^ Mq. As we expect, massive haloes produce much more 
systems with great velocity widths, while the less massive haloes 
produce more smaller velocity systems. 



5.2 Correlation between Metallicity and Velocity Widths 



Observational data by [Ledoux et al.| |2006| and [Prochaska et al.| 
( |2008^ confirmed in the last years the existence of a positive cor- 
relation between low-ion velocity width and the metallicity of 
DLAs along a given sightline. We showed that our reference runs 
fail to reproduce the DLAs velocity width distribution so we use 
the two post-processing criteria presented in the previous Section 
(MDW_rHi and 8^20, 512 smoothed) to test metallicity-velocity 
widths correlation. In our simulations we compute the metallici- 
ties of DLAs by taking the average metallicity value for the given 
LOS in a region of ±100 km/s in redshift space, centered at the 
location of the CM of each halo hosting a DLA. 

In the left panel of Figure [TS] we plot metallicities vs. ve- 
locity widths for 100 haloes randomly selected in a subsample 
of haloes with mass M > lO^°/i"^M0 in the MDW -Thi (green 
solid squares) and SW 20, 512 smoothed (red solid diamonds), along 
with observational data by [Prochaska et aL] ( |2008[ > (black crosses). 
Green squares (MDW .thi) show a weak correlation while cor- 
relation is negligible for the red diamond {SW2o,5i2smoothed). 

20, 512 smoothed have in average both metallicities and veloc- 
ity widths higher than MDW .thi and this is not surprising due to 
the procedure of smoothing metallicity over a large scale. For this 
plot we use only haloes with masses above M > lO"'/i"^M0 
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Figure 18. Left Panel: relationship between metallicity and Sill velocity width of individual sightlines thi'ough haloes in MDW.rjji (green solid squares) 
and SW 20, 512 smoothed (red solid diamonds) at 2 = 3, compared to the observational data by Prochaska et al. (20081 (black crosses). Right Panel: 
relationship between metallicity and Sill velocity width of individual sightlines through haloes with masses M > lO'-'-''^h~^ Mq (red solid diamonds), 
M > 10^^h~^ Mq (green solid squares) and M > W^'^h~^ Mq (blue solid triangles) in MDW_thi run at z = 3, compared to the observational data by 
[Prochaska et al.|f2008) (black crosses). 



because we noted that for the most massive haloes the correla- 
tion is stronger. In the right panel of Figure [TS] we clarify this 
point using the MDW.thi run. Red diamonds refers to all the 
haloes with mass Al > 10^'''^ Mq, green squares are the same 
green points of the left panel (a subsample of haloes with mass 
M > IO^^/i^^Mq), while blue triangles derive from a subsample 
of haloes with mass Af > lO^'^/i~'^M0. It is evident that the corre- 
lation in this case is somewhat weaker. Instead, using a subsample 
of haloes at larger masses, we increase the value of Avgo, since 
faster winds are produced by more massive haloes, that usually are 
more metal rich at least when compared with smaller mass haloes. 
Thereby, if we restrict ourselves to haloes above lO^'''^/i~^M0the 
agreement is reasonably good, even though smaller mass haloes, as 
we saw in the previous sections, contribute to all of the other DLA 
statistical properties. We also note that in general the mean metal- 
licities recovered from our simulations are somewhat smaller than 
the observed ones. 



6 DISCUSSION AND CONCLUSIONS 

In this paper we investigated the properties of DLAs in high- 
resolution hydrodynamical cosmological simulations. In particu- 
lar, we focussed on the role of feedback in the form of galactic 
winds to quantify their impact on the neutral hydrogen and metal 
ion species content both on the global IGM and on the IGM around 
putative DLA sites. The analysis was made by comparing with the 
recent works of Nagamin e et al!] ( |2007| > and |Pontzen et al.|p008| l, 
but with important differences with respect to these previous analy- 
ses: i) higher resolution simulations (by a factor ~ 10 in mass) than 
Nagamine et al. ( 2004[|2007 1; ii) statistics performed in a cosmo- 
logical setting, like in|Nagamine et aL] ( [2004[[2007^ , but differently 
from |Pontzen et al.]j2008^ , who relied mainly on simulations of sin- 



gle objects; in) inclusion of metal cooling and an accurate chemo- 
dynamical code that follows self consistently the IGM enrichment 
(Tomatore et al. 2007 >; iv) inclusion of different wind prescriptions 
(energy driven and momentum driven); v) modification of the stel- 
lar IMFs and the linear dark matter power spectrum in the initial 
conditions to explore the parameter space further. 

The conclusions that we draw from our analysis can be sum- 
marized as follows. 

- The different feedback prescriptions that we explored (WW, 
weak energy-driven winds of 100 km/s; SW, strong energy-driven 
winds of 600 km/s; MDW, momentum driven winds) give distinct 
predictions for the gas distribution in metallicity-temperature and 
density. SW are effective in heating the gas particles at tempera- 
tures of 10^ K to a larger extent than MDW and WW. Also the 
metallicity-density distribution is different: WW and MDW show 
correlations in high-density regions (denser regions being more 
metal rich), while SW is more efficient in polluting the low-density 
IGM (Figures 1 and 2). 

- WW, SW and MDW have also an impact on global quantities 
of the simulations such as the cosmic star formation rates and the 
neutral hydrogen content. As for the latter, MDW and SW show 
agreement with data while in the WW the HI is clearly overpro- 
duced since the feedback is not as efficient. Also, the star formation 
rate for the WW run is higher than for the SW one (see Figures 3 
and 14). 

- The evolution with redshift of CIV and OVI (two of the most 
common ion species observed in QSO absorption spectra) is similar 
in shape between the different runs but the normalization differ by 
up to a factor 3 (see Figure 4). 

- Focussing on the properties of the haloes that could host DLAs 
we find that haloes between 10^ — 10^'^h~^ Mq at 2: = 3 have sim- 
ilar SFR for the different wind models of about 0.01 — 0.1 Afo/yr, 
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with large scatter, while more massive haloes have larger SFRs and 
the trend with mass depends on the specific wind model (Figure 7). 
The metallicities of the haloes as a function of mass is correlated 
for haloes of masses > IO^'^/i'^Mq in the WW and MDW runs, 
while no correlation or possibly a weak anti-correlation is seen for 
the SW model. For less massive haloes below lO^'^/i"^ Mq the 
scatter in metallicity is huge and span several orders of magnitude 
(Figure 8). 

- The cross-sections inferred from the different runs are in over- 
all agreement with the results of |Nagamine et al.| ( |2004[ |2007| l, 
although the values for shape and normalization are somewhat 
smaller for similar runs. This is likely to be due to the different 
prescriptions of the wind implementation and the metal cooling. 
The DLA incidence rate is in good agreement with the observed 
one by |Prochaska et al.| ( (2005l l from SDSS data for MDW and SW, 
while is larger for WW. The incidence rate results clearly show that 
all the haloes whose masses are above IVph^^ Mq contribute to 
the cross-section (Figures 10 and 11). 

- The column density distribution function at z = 3 is in rough 
agreement with the data points for all the models (for MDW the 
agreement is better), while at lower redshifts there are larger differ- 
ences and SW (WW) underpredicts (overpredicts) the number of 
DLAs. The WW run shows the largest discrepancies especially for 
the high column density DLAs, a feature that was already found 
by lNagamine et ar| (2004 l. The contribution of haloes of masses 
between 10^ — 10^'^/i~^Mq is particularly relevant for DLAs be- 
low A*'hi ~ lO^'''* cm~'^, another hint that less massive haloes 
contribute to reproduce DLAs statistical properties (Figures 12 and 
13). 

- Qualitatively, physical quantities interpolated along LOS that 
pierce the haloes show different behaviours (particularly in temper- 
ature and metallicity) in different models. Also the peculiar velocity 
fields look different (Figure 15) and the gas velocity gradient is in 
agreement with the observed wind velocities. 

- The distribution of velocity widhts of low ioniziation species 
(in our case Sill) reproduces the observed one l |Prochaska et al.| 
|2008 [) only for velocity widths below 100 km/s, while larger values 
are dramatically under-reproduced by our simulations. A slightly 
better agreement is obtained when Sill is supposed to faithfully 
trace the HI distribution (as done by jPontzen et al.|2"008^ , but even 
in this case the observed distribution is underestimated. We can fit 
the observed velocity widths well only if we empirically smooth 
the metallicity of each gas particle over a region of 500 comoving 
/i~^kpc (Figure 16). We also find a trend of Aiigo with the halo 
mass that goes in the expected direction: more massive haloes pro- 
duce the larger values for Augo (Figure 17). 

- The metallicity-velocity width correlations is broadly repro- 
duced by the two models discussed above, even if the metallic- 
ity values are somewhat lower than the observed ones. If we split 
the contribution to the correlation by taking into account different 
ranges of masses, we find that there is good agreement if we restrict 
ourselves to the most massive haloes with M > W°-^h~^MQ 
(Figure 18). This is in agreement with |Bames & Haehnelt|p008| l, 
but in disagreement with other statistical properties of DLAs, that 
in order to be reproduced need less massive haloes as well. 

Overall, we succeeded in reproducing most of the observed 
properties of DLAs in particular the column density distribution 
function and the incidence rate. It appears that the best agreement 
is given either by the SW or the MDW implementation, while the 
WW does not seem to fit the data. For the metal distributions, as 
traced by the velocity widths, the agreement is not good. This could 



be the hint that there are some pieces of physics missing in our 
treatment. Possible ingredients that are not considered here and that 
could help in easing the discrepancy are: i) radiative transfer effects 
in DLAs that have been either neglected so far (e. g.^Nagamine et aT] 



|2004[ l or approximately modelled (e.g. [Pontzen et al.|20()8 1: even 
if their impact seems to be not dramatic at least using our effective 
model description of the ISM, since the cold clouds are assumed 
to be fully self-shielded and the ambient medium to be optically 
thin, they could be important for large column density systems, es- 
pecially when affecting the molecular hydrogen in star forming re- 



gions I Gnedin et al. 2008 i; ii) small-scale turbulence and its impact 



on the metal diffusion at large scales, which has been recently in- 
vestigated by |Scannapieco & Briiggen ( 2008 1 and that can be effec- 
tive in smoothing the metal distribution around DLAs; Hi) a pre- 
enrichment of the IGM, possibly produced by PopIII stars, at higher 
redshift than those considered here that could result in a smoother 
metallicity for the particles that are in a wind phase (e.g. |Tornatore| 
et al. 2007). We stress that all these effects are not considered in 
our simulations and their accurate description could be important 
to match the remaining still unexplained DLA properties. 
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